Skip to content

Draft UWTN 2026-016: Boundary conditions on non-planar boundaries - #26

Open
lmoresi wants to merge 30 commits into
mainfrom
note/non-planar-boundary-conditions
Open

Draft UWTN 2026-016: Boundary conditions on non-planar boundaries#26
lmoresi wants to merge 30 commits into
mainfrom
note/non-planar-boundary-conditions

Conversation

@lmoresi

@lmoresi lmoresi commented Aug 18, 2026

Copy link
Copy Markdown
Contributor

R1 in WRITING-PLAN.md, rescoped from "rotated boundary conditions" — the question is the better frame and it is the one a reader arrives with. At review.

"No flow through this wall" is a velocity component on a box, and is not a component of anything on a sphere, a deformed mesh, or a surface with topography.

What it says

Three approaches, in the order in which each answers the previous one's weakness, with the weak forms written out:

enforces leak
direct penalty weakly; consistent only in the limit set by the parameter — and closing it conditions the operator worse
Nitsche weakly; consistent above a stability threshold still a discretisation-order leak
rotated DOFs strongly machine precision

The maths is there because the differences are the argument — "add two more terms and it becomes consistent" is not believable unseen.

The rotated approach is presented as the classical answer it is. Engelman, Sani & Gresho were reviewing the alternatives in 1982, so the question the note actually answers is why the exact method is the least used of the three. That answer is structural rather than numerical: a per-node rotation leaves the problem in a mixed basis, and everything downstream has to agree about which nodes are which.

The figure shows where that obligation falls, and the point is that it is contained — the velocity solve is rotated and carries its multigrid with it, while the Schur complement and pressure solve abutting it never handle a rotated vector. One un-rotation on the boundary between them feeds both.

Also: a coordinate-system formulation is treated as the same rotation applied globally rather than as a fourth approach, which explains both its appeal and why it stops working the moment the boundary leaves a coordinate surface.

Two things left deliberately visible

  • The "when the choice matters" section carries a warning that it is not measured. The claim — that a traction recovered from an approximately satisfied constraint inherits the approximation, while the rotated reaction is σ_nn — is structural and we believe it. That is not the same as having shown it, and the note says so rather than asserting a number. This is the main outstanding work: one experiment, the same model under all three treatments, surface traction against a known answer.
  • The note asks whether Behr (2004) reaches the same measure-weighted normal we did. Its bibliographic details are verified; it has not been read. If it does, underworld3#560 was a rediscovery and the note should credit it as the source. Dropped as it is basically speculative as written (literally - "not read this, but maybe")

Notes

  • The figure migrates from underworld3, where it was parked because this article did not exist. Its .typ, generator and JSON geometry ship in examples/, with build instructions in the source header.
  • SVG normalised through scripts/rebuild_figures.normalise_svg — Typst emits 652 random glyph ids, which test-unit rejects.
  • All three DOI citations verified to resolve to the right papers before use. The first Engelman DOI I tried resolved to a different paper in the same volume.
  • Carries the allocator fix from UWTN 2026-015: Testing a solver against exact solutions #22. Without it new_article.py issued 2026-012, which is already gthyagi's — the very collision that fix exists to prevent, because it is not yet merged.

validate clean, test-unit 130 passed, PDF builds to 5 pages.

Underworld development team with AI support from Claude Code

R1 in the writing plan. 'No flow through this wall' is a velocity component on
a box and is not a component of anything on a sphere, a deformed mesh, or a
surface with topography, and the note is about what to do instead.

Three approaches in the order in which each answers the previous one's
weakness: a direct penalty, which is consistent only in the limit so closing
the leak means conditioning the operator worse; Nitsche, which carries the
boundary traction terms and is consistent above a stability threshold; and
rotating the degrees of freedom, which imposes the constraint exactly. The
weak forms are written out, because the differences between them are the
argument and 'add two more terms and it becomes consistent' is not believable
unseen.

The rotated approach is presented as the classical answer it is -- Engelman,
Sani and Gresho were reviewing the alternatives in 1982 -- so the question the
note actually answers is why the exact method is the least used of the three.
The answer is structural: a per-node rotation leaves the problem in a mixed
basis, and the figure shows where that obligation falls. It is contained. The
velocity solve is rotated and carries its multigrid; the Schur complement and
pressure solve beside it never handle a rotated vector.

The figure moves here from the underworld3 repository, where it was parked
because this article did not exist yet. Its typst source, generator and JSON
geometry ship in examples/.

Two things are deliberately visible rather than quietly resolved. The section
on when the choice matters carries a warning that it is not yet measured -- the
argument is structural and we have not run the experiment. And the note asks
whether Behr (2004) reaches the same measure-weighted normal we did, because if
it does, our #560 was a rediscovery and should say so.

Carries the allocator fix from the exact-solutions branch. Without it
new_article.py issued 2026-012, which is already gthyagi's -- the same
collision the fix exists to prevent, because the fix is not yet merged.

Underworld development team with AI support from Claude Code
@github-actions

github-actions Bot commented Aug 18, 2026

Copy link
Copy Markdown
Contributor

Preview

Or the whole site.

Built from 1cc96b4. Shows notes at draft and review, which the published site withholds. Not indexed, no comments, and not the citable version. Confirmed serving this commit before this was posted.

The note asserted that a weakly imposed constraint holds to the accuracy of the
discretisation and a strong one to the accuracy of the arithmetic. Annulus, no
slip inside, the treatment under test outside, degree-four radial forcing;
largest normal velocity on the outer boundary against the true radial
direction, over the flow speed:

  cell    Nitsche    rotated
  0.150   4.6e-03    7.3e-11
  0.100   2.2e-03    6.5e-11
  0.075   1.7e-03    1.0e-10
  0.050   5.8e-04    8.4e-11

Refinement is what separates them and one resolution would not have. Nitsche
falls as about h^1.9 -- discretisation-limited, which is what consistency buys
and all it buys. The rotated constraint does not move with the mesh at all.

Three things the measurement cost, all worth recording:

The first control was worthless. Two free-slip circles leave the rigid rotation
unconstrained, and that nullspace is purely TANGENTIAL, so a radial leak metric
read 2e-14 on a solve that had diverged with |u| = 2.7e5. The inner boundary is
now no-slip and the control -- outer boundary left natural -- reads 0.98.

A vector MeshVariable's .array is (N, 1, dim). It broadcasts silently against
an (N, dim) array of normals and returns projections of 1e-16 for a velocity of
1e-2, without raising. Every number in the first run was that artefact.

Underworld exposes no separate direct-penalty condition and building one from a
natural BC would not converge, so that row has no measured counterpart. Said in
the note rather than quietly dropped, alongside the other limit: this measures
the constraint, not the traction recovered from it.

Underworld development team with AI support from Claude Code
The previous commit said Underworld exposes no direct-penalty condition. It
does, it is documented in docs/advanced/curved-boundary-conditions.md, and it
is in a good many working scripts:

    G = mesh.Gamma
    stokes.add_natural_bc(1.0e4 * G.dot(v.sym) * G, "Upper")

What I had written used the nodal boundary_normal and a NEGATIVE coefficient.
A negative one is anti-damping, so the linear solve failed on the first
iteration, and I read that as the method being unavailable rather than as my
sign being wrong.

With all three measured, the resolution table gains its first column, and the
two weak methods look alike in it -- both leak parts in a thousand and both
improve with the mesh. What separates them is their own parameter:

  penalty  1e2 2.6e-01 | 1e3 2.3e-02 | 1e4 8.5e-04 | 1e5 9.6e-04 | 1e6 diverged
  gamma      1 diverged |  10 1.7e-03 |  100 2.7e-04 | 1000 3.0e-05

The penalty improves in proportion to how hard it pushes until the conditioning
catches up: no gain from 1e4 to 1e5, failure at 1e6. Nitsche fails below its
stability threshold and improves steadily above it, with no wall in this range,
which is why gamma=10 can be a documented default instead of a number to fit
per problem.

Two runs in the first parameter sweep had failed the line search and were about
to be tabulated as measurements. The sweep now records the converged reason and
prints 'diverged' instead of a number.

Underworld development team with AI support from Claude Code
The note said Nitsche 'keeps improving without the conditioning wall'. That was
written from a sweep that stopped at gamma = 1000, which is one value short of
where it fails.

  gamma    1e3  3.0e-05  converged
           1e4  --       DIVERGED_LINE_SEARCH
           1e5  --       DIVERGED_LINE_SEARCH

It does not stop converging in the metric, it stops converging: reason -6, and
|u| collapses toward zero, so the leak figures above the threshold are from
solves that never finished.

So both weak methods are bounded above by conditioning and neither escapes
tuning. What differs is the floor each reaches before failing -- the penalty
bottoms out near 1e-3, Nitsche at 3e-5 -- which is what consistency actually
buys: a lower floor, not an unbounded parameter.

Nitsche is bounded below as well, at gamma = 1, where the form stops being
coercive. Its usable range has a threshold at each end, and gamma = 10 sits in
the middle of that window on any mesh because gamma is dimensionless and the
term it scales already carries mu/h. The penalty coefficient carries no such
scaling, which is why its working value is a property of the problem.

The shipped sweep now runs into the failing range rather than stopping just
below it.

Underworld development team with AI support from Claude Code
The note had three approaches and there are four. add_constraint_bc carries a
scalar multiplier field h coupled into the saddle point, so the constraint is a
ROW of the system rather than a term added to one -- and at convergence h on
the boundary is the normal traction, which makes it the second of the four that
returns the stress rather than recovering it.

Measured alongside the others, same annulus, same forcing:

  cell    penalty   nitsche   multiplier  rotated
  0.150   4.5e-03   4.6e-03   8.3e-04     7.3e-11
  0.100   2.5e-03   2.2e-03   1.6e-04     6.5e-11
  0.075   8.5e-04   1.7e-03   5.6e-05     1.0e-10
  0.050   9.5e-04   5.8e-04   1.2e-05     8.4e-11

It is an approximation, as expected, but a much better behaved one: it starts an
order of magnitude below the weak forms and falls near h^3.9 against their
h^1.9, because the only thing approximated is the screening of the interior
multiplier degrees of freedom rather than the enforcement itself. The note says
where that screening enters, and that the augmented-Lagrangian term conditions
the [p,h] block without biasing the multiplier, so accuracy does not depend on
it the way a penalty's does.

Also started the surface-traction comparison, which is the test that matters,
and it is NOT in the note. Two exact routes to sigma_nn -- the multiplier and
the rotated reaction -- agree to about 1% at one resolution, but under
refinement the vertex agreement converges while the midpoint agreement gets
steadily worse, which is backwards: the solver documents midpoints as
superconvergent on a curved boundary. The fault is most likely this script
selecting the multiplier's boundary trace by a radius band that narrows with
the mesh; it needs to come from the boundary label instead. The symptom and the
suspicion are recorded in examples/traction.py rather than the numbers being
quoted anywhere.

Underworld development team with AI support from Claude Code
The traction comparison in traction.py has no oracle -- it measures two methods
against each other, so its diverging midpoint column cannot be interpreted.
uw.analytic.Zhong2008 fixes that: .response() returns the surface topography
kernel directly, 0.4191904 for the default degree-2 case, referenced to Zhong
et al. (2008) GGG 9 Q10017 after Hager and O'Connell. The quantity it publishes
is the one the note cares about rather than a velocity field to infer a stress
from.

One consequence worth stating before anyone starts: Zhong is SPHERICAL. The
harmonic degree is a spherical-harmonic degree and the solution carries a planet
radius, so the matching model is a 3-D shell rather than the 2-D annulus the
leak tables use. That sets the cost and is why this is separate work rather than
another column.

Also records why SolCx stays useful even though it cannot test a curved
boundary: on a box all four treatments reduce to the same component constraint,
but it is still a fair test of whether a penalty or Nitsche form constrains the
normal degrees of freedom it is handed, with the geometry factored out and an
exact answer to hand.

Underworld development team with AI support from Claude Code
Checked what the analytic suites actually vary. assess's cylindrical and
spherical solutions (Kramer et al. 2021) take nu as a SCALAR -- isoviscous.
Zhong2008 layers viscosity by radius through viscosity_interfaces, also not
laterally. SolCx has the lateral jump and a Cartesian box.

So no analytic solution has both a curved boundary and a lateral viscosity
variation, and the test splits in two rather than one being a weaker version of
the other: Kramer or Zhong for the geometry with the rheology trivial, SolCx
for the rheology with the geometry trivial. Lateral contrast is the condition
under which these methods are known to give trouble and it is the only one of
the two for which an oracle exists at all.

Also noted that a shell with laterally varying viscosity has no exact solution,
so that case needs a different kind of evidence -- convergence or agreement
between methods -- and the note should say which rather than implying a
measurement it cannot make.

Underworld development team with AI support from Claude Code
The leak tables measure the constraint. This measures what the constraint is
wanted for, against `uw.analytic.CylindricalStokes` (Kramer et al. 2021), whose
`assess` package publishes the radial stress rather than leaving it to be
recovered. On the outer boundary the exact answer is a single harmonic,
0.1506696 cos 2t, fitted to a residual of 1e-16, so the metric is one amplitude.

Three findings, two of which correct what the note said.

The direct penalty written against the FACET normal does not converge. Its leak
falls to 1e-5 while the velocity is 60% wrong and the surface stress 26% wrong,
and refining from cell 0.15 to 0.035 does not move either. The same coefficient
against the measure-weighted node normal converges cleanly. This is the classic
over-constraint the consistent normal was introduced to avoid, and it means the
note's penalty column was measuring the onset of locking: at a coefficient of
1e3 that penalty leaks 3e-2 and gets the stress right to 2e-3, and at 1e8 it
leaks 1e-7 and is 26% wrong. Both penalty columns are now tabulated.

Once the constraint is imposed against the node normal, which treatment imposed
it stops mattering to the recovered stress -- rotated, multiplier, Nitsche and
penalty all read 6.3e-3 at cell 0.075. The recovery sets the floor, not the
boundary condition, so the note's reasoning that a weakly imposed constraint
poisons the traction was too strong. What the rotated constraint buys is the
reaction, which is three times better on the same solve and differentiates
nothing.

The multiplier is not quite the whole traction: the momentum row carries
h + r(u.n - g), and with the default r = 1e4.mu the missing share is a few per
cent of the surface stress. underworld3#607.

Also measured: the SolCx half, where Nitsche at gamma = 10 gets the surface
stress 277% wrong at a viscosity contrast of 1e6 while holding the constraint to
1e-7, and where the rotated constraint disagrees with a component condition at
the corner they share (underworld3#608).

Behr (2004) is read and cited for what it says: the consistent normal is
standard and preferred for conservation, and it does not on its own guarantee a
good discrete slip condition. The derivation in the "Which normal" section is
credited to Engelman, Sani and Gresho (1982) rather than presented as ours.

`traction.py` is removed -- `stress.py` supersedes it, with an oracle and with
the boundary trace selected by label rather than by a radius band. The mesh
cache the examples write beside themselves is untracked and gitignored.

Underworld development team with AI support from Claude Code
@lmoresi

lmoresi commented Aug 19, 2026

Copy link
Copy Markdown
Contributor Author

The stress test is in. It measures the surface normal stress against an exact
answer rather than the constraint, and it changes two things from the original text.

The facet-normal penalty does not converge. Its leak falls to 1e-5 while the
velocity is 60% wrong and the surface stress 26% wrong, and refining from cell
0.15 to 0.035 leaves both where they are. The same coefficient against the
measure-weighted node normal converges cleanly — one line of difference. So the
note's original penalty column was measuring the onset of locking rather than
compliance: at a coefficient of 1e3 that penalty leaks 3e-2 and gets the stress
right to 2e-3; at 1e8 it leaks 1e-7 and is 26% wrong. Both penalty columns are
now in the tables, and there is a figure.

The recovery is the floor, not the boundary condition. At cell 0.075 the
rotated constraint, the multiplier, Nitsche at gamma=100 and the penalty at 1e5
all recover the surface stress to 6.3e-3. The argument that a weakly imposed
constraint poisons the recovered traction is not what the measurement shows once
the constraint is written against the right normal. What the rotated constraint
buys is the reaction — 2.1e-3 on the same solve, and nothing is differentiated
to get it.

Two defects surfaced, both filed against underworld3:

  • #607 — the
    multiplier is not the whole traction. The momentum row carries
    h + r(u.n - g); with the default r = 1e4·mu the missing share is a few per
    cent of the surface stress, and the docstring says accuracy is independent
    of r.
  • #608 — a rotated
    free-slip wall and a component Dirichlet wall disagree at the corner they
    share, by 0.23 in the pressure at that one node, with everything else agreeing
    to 4e-4.

The SolCx half is in too. It cannot rank the treatments — the component
Dirichlet control reads 0.085 with a velocity error of 8.8e-6, because on the
stiff half the recovered stress is a cancellation against 2·1e6·du/dz — but it
does say that Nitsche's gamma has to be raised with the viscosity contrast: at
1e6 contrast, gamma=10 gets the surface stress 277% wrong while holding the
constraint to 1e-7, and gamma=1000 returns it to the control floor.

Behr (2004) is now read and cited for what it says — the consistent normal is
standard and preferred for conservation, and does not on its own guarantee a
good discrete slip condition — and the "Which normal" derivation is credited to
Engelman, Sani and Gresho (1982) rather than presented as ours.

Also: a banner, traction.py removed in favour of stress.py (oracle, and the
boundary trace selected by label), and the examples' mesh cache untracked.

Still open: the spherical case (Zhong 2008, surface_topography = 0.41919),
which is the 3-D cost step.

Underworld development team with AI support from Claude Code

lmoresi added 19 commits August 19, 2026 17:39
Predicted against computed, along the top wall, at a viscosity contrast of 100
and of a million. The exact answer is published by `uw.analytic.SolCx`, so this
comparison can be a figure rather than a column of numbers.

Drawing it turned up a defect the tables had missed. At a contrast of 1e6 the
multiplier field as `multiplier()` returns it peaks at 0.042 against an exact
0.383 and is anti-correlated with it (-0.53), on a solve whose velocity error is
8.8e-06. The default augmentation is r = 1e4.mu(x), which the viscosity step
makes 1e10 on the stiff half, so almost the whole traction sits in the
augmented-Lagrangian term rather than in h. Adding it back gives 0.047, the most
accurate surface stress in the comparison -- better than the projected recovery
at 0.075. Turning r down is not available: at this contrast r = 0 and r = 1e2
both fail to solve. Recorded on underworld3#607, which the annulus had already
found at the few-per-cent level.

The figure carries both multiplier curves, as returned and corrected, so the
size of the missing term is visible rather than asserted.

The script caches its curves beside itself: the solves take about ten minutes
and the figure is redrawn more often than it is recomputed. A dash pattern is a
linestyle and not a format string, and a sign convention is asserted rather than
assumed -- the first draft of this figure stopped on that assertion, which is
how the collapse was noticed.

Underworld development team with AI support from Claude Code
… are masked

Two corrections and one identity, all from Louis reading the draft against the
free-surface work.

The note said the multiplier's interior degrees of freedom are screened and that
the screening is where the approximation enters. They are not screened, they are
CONSTRAINED OUT of the global system in the PetscSection before the solve sees
them (_constrain_interior_multipliers_in_section), which is documented as
lossless: the interior rows are the screening block alone, so the interior
multiplier is determined by the boundary trace and pinning it moves a converged
solve by ~1e-8 in velocity. They neither cost anything nor limit the accuracy.

The identity: the momentum row's boundary load is M_Gamma (h + r(u.n - g)), and
at convergence that balances the volume residual restricted to the boundary,
which is the nodal load the consistent boundary flux back-calculation reads
(Zhong, Gurnis & Hulbert 1993). So the corrected multiplier is the CBF traction
de-smeared by the boundary mass -- the same computation, arrived at by carrying
the traction as an unknown rather than reading it out of the residual. Measured
across two solves they agree to 3.2% and 4.9%, inside each route's own error.
That is why the free-surface work was right on both counts: it rejected the
multiplier AS RETURNED, which is missing the augmentation share, and kept the
rotated lid with the CBF reaction, which is the same quantity with nothing
missing.

The corner spike is re-framed. It is the documented behaviour of the reaction
recovery at a node where a rotated wall meets an essential one -- the value
reported there is the essential constraint's reaction, not the wall's traction --
and the practice, which the free-surface SolCx work already followed, is to
leave the corner out of the rotation. Control: the same CBF read on a lid held
by the component condition alone peaks at 0.381 against an exact 0.379 with no
spike, where the rotated run peaks at 0.497. underworld3#608 is re-scoped to
match.

New: underworld3#614, boundary_flux() returns ~1e12 on Stokes_Constrained while
reading the exact topography to 8% on an ordinary solve. It is the cross-check
that would verify the identity within a single solve, and it is unavailable.

The Zhong citation is pinned in a references.bib: the DOI contains parentheses
and MyST's inline DOI form truncates at the first one. An at-sign in a bibtex
COMMENT fails the whole file.

Underworld development team with AI support from Claude Code
…taminated table

The three defects this note turned up are fixed in underworld3 PR #617, so every
table here is re-measured against a solver that behaves.

The SolCx half is the one that changes. The three exact treatments -- component
Dirichlet, multiplier, rotated -- now agree to three figures at every viscosity
contrast from 10 to 1e6, whole wall and trimmed alike. Before, the rotated column
read 0.322 at a contrast of 10 against the reference's 0.048, all of it two
corner nodes, and the multiplier was an order of magnitude out at 1e6. Getting
those three to agree was the point.

A bare penalty coefficient now has numbers instead of "diverged": 0.045 at a
contrast of 10, and 0.992 at 1e6, which is to say none of the signal. Scaling it
by the local viscosity still does not solve.

WITHDRAWN: the earlier SolCx numbers, including the Nitsche column. They were
taken while concurrent runs shared one mesh-cache file -- StructuredQuadBox keys
its cache on the box corners and NOT on the element resolution, so runs at
different resolutions silently swap meshes (underworld3#618, filed). A marginal
solve then flips between converged and diverged for reasons that look like the
method, which is exactly how it fooled me. Everything in that section was
re-measured sequentially on a cleared cache. Nitsche has no column now: our
configuration of it on this box converges at a contrast of 1e6 and fails the line
search at 10, at both resolutions and at gamma = 10, 100 and 1000, and we are not
confident enough in that setup to publish numbers from it.

New: what each treatment COSTS. The solve is the same for all four within the
spread; the recovery is 15 to 25 times cheaper for the reaction and the
multiplier than for a projection, because they read the state the solve already
returned. And a "which one to use" section, which the accuracy and the cost now
point at together.

The multiplier section is rewritten around the fix rather than the defect, with
the clean SolCx numbers: h alone reads 0.042 against an exact 0.381 and is
anti-correlated with it; h + r(u.n - g) reads 0.382.

Underworld development team with AI support from Claude Code
The preview build failed on CI with "Citation data from doi.org was not available
or malformed" for the Behr and Kramer DOIs, while resolving fine locally. That is
the failure the MMPDE note's references.bib was written to prevent: a build that
depends on doi.org answering will eventually publish a note with a broken
citation, and a deposited PDF cannot be repaired afterwards.

All five are now pinned and cited by key.

Underworld development team with AI support from Claude Code
Louis's note: the interesting thing here is two penalty approaches against two
ways of imposing the constraint exactly -- by construction, or by a Lagrange
multiplier -- and that pairing is now stated where the four are introduced.

Removed: the section that told the multiplier's missing augmentation share as a
defect story, the withdrawn-measurement note, the "how it survived validation"
paragraph, and the issue numbers. What survives is the method. The traction a
multiplier constraint holds the boundary with is h + r(u.n - g), which is stated
where the method is described, and the two exact routes are shown to be the same
computation, which is the point worth keeping.

The timings are re-measured at sizes where they mean something: 71 000 velocity
nodes rather than 2 000, where the numbers were hundredths of a second and inside
the run-to-run spread. At that size the multiplier's solve costs 50% more than
the others, and the recovery separates by a factor of 25 to 40 -- the projection
costs MORE than the Stokes solve (0.63 s against 0.45 s) while the two exact
routes read their traction off the state the solve already returned (27 ms and
16 ms). The smaller meshes are kept only to show the scaling, and the table says
so.

Underworld development team with AI support from Claude Code
…er implies one

Louis asked why determining topography from the multiplier requires a solve. It
does not, and neither does the rotated reaction: the multiplier's traction is
h + r(u.n - g) evaluated on the boundary trace, and the rotated constraint's is
its nodal reaction divided by the LUMPED boundary mass, which is a division. In
2-D neither route assembles anything.

The 27 ms in the previous table was the example's bookkeeping, not the method.
Profiled at 71 424 velocity nodes with 1 008 trace nodes: 2.7 ms to walk the DM
for the boundary node list, 5.3 ms for the traction arithmetic, and 18 ms
building a KD-tree over the WHOLE mesh -- twice, once for the multiplier and once
for the velocity -- to line the trace up with the field arrays. The multiplier is
carried at the velocity degree, so one lookup serves both, and it depends on the
mesh rather than the solution: a time-stepping consumer builds it once. With that
fixed the recovery reads 9 ms against the rotated reaction's 16 ms.

So the table now says the thing that is actually structural: the weak forms have
to RECOVER their traction with a projection solve, which costs more than the
Stokes solve itself (0.63 s against 0.45 s), while the two exact routes read
theirs off the state the solve already returned in 9 to 16 ms.

Underworld development team with AI support from Claude Code
The multiplier has nothing to recover: h is a finite element field, so its nodal
values ARE the traction pointwise, and h + r(u.n - g) is an expression evaluated
where it is wanted. A millisecond figure for "reading a field" invited exactly
the comparison it should not -- the slot is a dash.

The rotated reaction does need a step, and the note now says which. The reaction
is an INTEGRATED nodal load, M_Gamma times the pointwise traction, so recovering
a pointwise value means undoing the boundary mass. On a 2-D trace and on 3-D P1
triangles that mass is lumped and undoing it is a division -- the 10 to 16 ms
measured here. On 3-D P2 triangles it is a genuine solve: the lumped row sums
vanish at the vertices, so the consistent trace mass is assembled and solved, and
Underworld gathers the trace to one rank to do it. That is the case a spherical
free surface runs in, so it is worth naming.

Which makes the real distinction not speed but what each method hands you: the
multiplier gives the traction as a field, the reaction gives it as a load that
still has to be divided by a mass.

Underworld development team with AI support from Claude Code
… method's

Louis asked whether the weak forms really need a global re-solve to get stresses.
They do not need a GLOBAL one -- an L2 projection over every node to get values
on a thousand boundary nodes is more work than the job requires -- but the cheap
route is not open to them, and measuring it says why.

The consistent boundary flux reads the assembled residual instead of
differentiating anything. On the rotated constraint it matches the projection at
a third of the cost (error 1.6e-4 either way). On Nitsche and on the penalty it
returns an error of 1.00, which is the metric reporting that nothing was
recovered.

A reaction exists in the residual only where a row has been CONSTRAINED. A weak
condition supplies its traction as a term inside the row it acts on, so the
residual there is balanced at convergence and there is nothing left to read. The
multiplier does exactly the same thing and gets away with it because the term it
supplies IS the traction as a field; Nitsche's term is written in sigma(u), so
reading it back still means differentiating the answer.

Which lines up with the two pairs: imposed weakly, the traction is a by-product
and you differentiate for it; imposed by construction, it is the constraint
reaction and you de-smear a nodal load; imposed by a multiplier, it is an unknown
of the system and you read the field. The cost of the first is negotiable, the
differentiation is not.

Underworld development team with AI support from Claude Code
The note opens with a free surface and now closes there. The point of the
recovery taxonomy, for that case, is that sigma_nn is needed EVERY step -- the
surface moves under the traction it carries -- so a recovery that costs more than
the Stokes solve doubles the model rather than adding a post-processing step at
the end of a run.

Stated for the two exact treatments, since they are what a free surface should be
choosing between. The rotated constraint's reaction on the 3-D P2 triangular
trace a deforming spherical surface actually uses is de-smeared by a
consistent-mass solve that is gathered to one rank: a serial section inside a
parallel timestep. The multiplier hands over a field, so there is nothing to
de-smear in 2-D or 3-D, and pays about 50% on the solve instead.

Which of those wins on a large parallel shell is a question of rank count and we
have not measured it. The note says so and names it as the next measurement
rather than guessing.

Underworld development team with AI support from Claude Code
Four things Louis asked for, on top of his edit.

The penalty term no longer carries gamma/h. It is a bare kappa, and the text says
so: none of the mu/h that Nitsche's parameter carries, which is exactly why the
value that works is a property of the model. The two statements contradicted each
other before.

The multiplier is lambda, not h. It was colliding with the mesh size in gamma/h
and with the topography it produces. The prescribed wall-normal velocity is
u-tilde-n rather than g, which was colliding with gravity. The section now says
what lambda IS -- it has units of stress, at convergence it is sigma_nn, and
dividing by delta-rho g (gravity, no relation to the constraint datum) gives a
dynamic topography.

The choice of normal is settled in its own section and not revisited. The
locking measurement and the figure move there, the separate "a constraint that is
satisfied, and wrong" section is gone, and the leak, parameter and stress tables
carry ONE penalty column -- the node normal, which is what everything else uses
by default. The section ends by saying the facet normal does not appear again,
and it does not.

The passage Louis marked unclear is replaced. With one normal in the results
there is nothing to disambiguate: the penalty trades a decade of coefficient for
a decade of leak, with no wall in this problem, and what it does not do is
converge.

Also fixed, from the edit: \label and \ref are LaTeX and KaTeX has neither, so the
build was failing 92 times on "Undefined control sequence: \ref". The equation is
labelled the MyST way and referenced with {eq}, which renders it as "(1)". That
needed one more thing to work -- a display-math block opening on the line
directly after prose is parsed as part of that paragraph, and its label is left
as literal text. Every $$ block in the article now has a blank line before it.

Underworld development team with AI support from Claude Code
The note was holding back the one thing a reader comes for. Each of the four
sections now ends with how you get the surface stress out of that method, with
the expression:

  penalty, Nitsche   sigma_nn = n.sigma(u,p).n recovered from the solution --
                     a projection, and a second solve -- because the boundary
                     term supplies the traction inside the row it sits in, so
                     the residual carries none
  multiplier         sigma_nn = lambda + r(u.n - u~_n), read on the trace: no
                     recovery, because lambda is a field
  rotated            sigma_nn = -M_Gamma^-1 (Au - b)|_Gamma, the reaction of the
                     struck row de-smeared by the boundary mass -- the CBF

and h = -(sigma_nn - mean)/(delta-rho g) in each case. The later section that
compared them no longer opens on "the correction above is not a patch", which
referred to something that had already been deleted; it opens by pointing at the
two expressions and showing they are the same one.

The benchmark is introduced once. It was described twice, in different words,
with the leak measured on one problem and the stress on another -- so the two
tables could not honestly be read against each other. They now come from THE SAME
SOLVES: `stress.py both` reports the leak and the stress error per run. The
numbers change slightly because the leak is now measured on the Kramer problem
rather than on a separate degree-four one.

Reading them side by side is the point, and it is a stronger result than either
alone: the rotated constraint holds the boundary nine orders of magnitude better
than the multiplier and lands on the same stress error, because what sets that
number is the recovery and not the boundary condition.

Citations are all one form now -- bracketed [@key], with the authors named in the
prose where the sentence wants them. One was still a raw DOI, which would not
have resolved against the pinned bibliography.

Underworld development team with AI support from Claude Code
The preview was failing while build and test passed, which is the signature of a
PDF-only fault: the HTML tolerates an unresolved reference and Typst does not.

    $$
    ...
    $$
    (eq-free-slip)          <- a paragraph of literal text

MyST reads the label only when it is on the same line as the closing $$. On its
own line the equation is never labelled, so the three {eq}`eq-free-slip`
references have no target: HTML renders them as bare spans and carries on, while
`typst compile` stops with "label <eq-free-slip> does not exist in the document"
and the archival PDF is never produced.

Verified by building the PDFs with UWTN_PREVIEW=1, which is what CI does -- the
note is at status: review, so a production build skips it and the fault does not
show. 45 PDFs, this one 766 kB, no label error.

Underworld development team with AI support from Claude Code
Louis's revision, brought in from his working copy. The substance of it: the
classical-strategy paragraph moves to where the rotated method is actually
described rather than trailing it, the normals section states what a discrete
normal IS before saying which one to use, the CITCOMS reference is given its
place in the coordinate-rotation discussion, and the forward references to
free-surface work that is not written yet are gone -- along with the passage
motivated by the traction requirement of our own free-surface algorithm, which a
reader cannot check and which belongs in its own note.

Two mechanical things carried with it.

The equation label needs a space: `$$ (eq-free-slip)`, not `$$(eq-free-slip)`.
Same failure mode as the label on its own line -- the equation goes unlabelled,
`{eq}` resolves to nothing, the HTML shrugs and `typst compile` stops.

Removing the Behr caveat left `Behr_2004` in references.bib with nothing citing
it, so the entry goes too. If that caveat should come back -- the consistent
normal not guaranteeing a good discrete slip condition, with recirculation at
curved walls -- the entry comes back with it.

Verified before pushing: 137 tests, a full UWTN_PREVIEW build with no MyST errors,
the archival PDF at 758 kB, and the equation rendering as a numbered, linked (1).

Underworld development team with AI support from Claude Code
Clarified the weak form of boundary conditions and introduced a direct penalty term for flow through boundaries.
Mostly just adding newlines around $$ so they render correctly. I hope I got them all !
The topography figure drew the multiplier's traction at twice its augmentation
share, so at a 1e6 contrast it left the panel -- exactly the "out by 2" it looked
like -- while the curve labelled "as returned" was the correct one.

The cause was two copies of one expression. `solcx.reaction_traction` was changed
to return `traction()` (lambda + r(u.n - u~_n), the whole boundary load, from the
solver) when underworld3#617 landed, and this script went on adding r(u.n) to
whatever that returned, which had been right while it returned the bare
multiplier. Nothing said the two had to agree, so they drifted the moment one
moved.

Both curves now come from the solver and neither is assembled here:
`multiplier()` for lambda, `traction()` for the whole load. The local copy of r
is gone, which is the part that could go stale. Recomputed from scratch: at 1e2
lambda 0.3822 and the traction 0.3811 against an exact 0.381; at 1e6 lambda 0.0415
and anti-correlated, the traction 0.3815 at corr +0.999. The labels say which is
which rather than "as returned", since what the API returns is now the traction.

The cetz figure says "A deformed surface has no preferred direction". Its SVG is
recompiled and run through `normalise_svg`, so Typst's per-compile random ids do
not dirty the tree.

Underworld development team with AI support from Claude Code
…ons' into note/non-planar-boundary-conditions
It was in the file and in nobody's way: I dropped the `{figure}` block when the
two benchmark introductions were merged into one, so `figures/topography.png` was
generated on every build and included by nothing. Which is why it was not in the
preview.

It belongs beside the contrast table in the SolCx section, and its alt text and
caption now describe the corrected figure -- lambda alone flat near zero at 1e6,
the traction on the exact curve, and the penalty failed by that contrast. The
figure count in metadata.yml goes to three.

Verified: full UWTN_PREVIEW build, no MyST errors, all three figures fingerprinted
into the page, archival PDF at 878 kB.

Underworld development team with AI support from Claude Code
Nitsche has no column in the SolCx table. Say what we measured -- our
implementation does not solve a boundary that mixes essential patches with
Nitsche patches, at any penalty we tried, while all four walls imposed weakly
converges -- instead of the parameter sweep and the speculation about which
way round the failure ought to have gone.

Confirmed at gamma = 100 as well as 10 and 1000: top wall only diverges in
the line search (SNES -6) at contrasts of 10 and 1e6, all four walls converges
(SNES 3), velocity error 1.74e-3 and 2.6e-5 against the exact solution.

The traction paragraph becomes a MyST note admonition, in the house form used
elsewhere in the corpus, with the claim as the title rather than as the first
words of the paragraph.

Underworld development team with AI support from Claude Code
Alexander Londono's aerial photograph of a skate bowl, cropped from the
3000x3854 original to 2400x840 -- the banner box's own 2.86:1, so the page
crops nothing further and the skater stays where the crop put him, on the
curved lip with his shadow beneath him.

Replaces the rendered annulus pair that stood in for it. The Unsplash credit
goes in metadata.yml, from where banner_body.py rebuilds it into the page
alongside the image.

Underworld development team with AI support from Claude Code
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

1 participant